Partie I - Initialisation:¶
Import des packages utilisés¶
import requests
import urllib.request
from bs4 import BeautifulSoup
import os
import pandas as pd
import numpy as np
from sklearn import linear_model
from sklearn import metrics
import matplotlib.pyplot as plt
import copy
from sklearn.linear_model import LassoCV
from netCDF4 import Dataset
from datetime import datetime
from datetime import timedelta
import plotly.express as px
import math
import warnings
warnings.filterwarnings('ignore')
Paramètres d'éxécution:¶
download_last_datas=False
#if selecting True the execution of some cell can be very long. In order to avoid it the results of thoose cells are
#stored in file and loaded if False is selected here
launch_long_calculations=False
#Lancement du dernier graphe dynamique, très long à executer
launch_last_graph=False
if launch_long_calculations:
from global_land_mask import globe
Définition des variables globales¶
global classement_url
global glossaire_url
classement_url="https://www.vendeeglobe.org/fr/classement"
glossaire_url="https://www.vendeeglobe.org/fr/glossaire"
Partie II - Scrapping et mise en forme des données¶
Fonctions de scrapping et de mise en forme du dataframe¶
def get_files_names(classement_url):
#list all the excel files of ranking in the classement page
r=requests.get(classement_url)
return [val.attrs.get("value") for val in BeautifulSoup(r.content.decode('utf-8'), 'html.parser').find("select",{"class":"form__input m--select onsubmit-rank"}).findAll("option") if val.attrs.get("value")!=""] if r.status_code==200 else []
def download_all_files(files_names):
#download all the files on the ranking pages. Delete the folder if exist and recreate it.
delete_if_exist_and_recreate_directory(os.getcwd()+"\\downloaded_rank_files\\")
for file in files_names:
url='https://www.vendeeglobe.org/download-race-data/vendeeglobe_'+file+'.xlsx'
urllib.request.urlretrieve(url, os.getcwd()+"\\downloaded_rank_files\\"+file+'.xlsx')
return files_names
def make_dict_boat_infos(boat):
# Function that create a dictionnaire over a boat (here a list getted with Beautifullsoup)
boat_dict={}
for item in boat:
boat_dict.update({item[0][:item[0].find(":")-1]: item[0][item[0].find(":")+2:]})
return boat_dict
def clean_boat_dict(list_boat_dict):
#Reagence the dictionnaire using all the boat dictionnaire created before, create a dictionnaire
#where the key is the number of the boat, and each value is the other element of the dictionnaire
#As the numéro de voile doesn't always match beetween the ranking files and the glossaire page, change the ones
#that are different in order to match correctly
final_dict={}
for boat in list_boat_dict:
nv=boat.get('Numéro de voile')
if nv!=None and nv!='ITA 06':
nv=nv.split(' ')[-1]
del boat['Numéro de voile']
try:
nv=int(nv)
if nv==16:
nv=10
except:
nv=int(nv[3:])
if nv==77:
nv=777
elif nv==None:
nv=59
else:
nv=34
final_dict.update({nv:boat})
return final_dict
def request_boats_info(glossaire_url):
#Function that use the different precedeing function in order to get the clean boat dictionnary
r=requests.get(glossaire_url)
return clean_boat_dict([make_dict_boat_infos([info.contents for info in boat.children if info!='\n']) for boat in BeautifulSoup(r.content.decode('utf-8'), 'html.parser').find("section",{"class":"boats-list"}).findAll("ul", {"class":"boats-list__popup-specs-list"})]) if r.status_code==200 else []
def delete_if_exist_and_recreate_directory(path):
#small function that check if a folder exist, then delete it.
#After create it again (or for first time)
try:
os.rmdir(path)
except:
pass
try:
os.mkdir(path)
except:
pass
return None
def download_files():
#Function that merge the two preceding function in order to download file more conveniently
return download_all_files(get_files_names(classement_url))
def get_df_from_excel(file, add_extension=False):
if add_extension:
path=os.getcwd()+"\\downloaded_rank_files\\"+file+'.xlsx'
else:
path=os.getcwd()+"\\downloaded_rank_files\\"+file
df=pd.read_excel(path, header=3)
first_line=df[df.index==0].values
df.columns=[first_line[0,i] if df.columns[i][:7]=="Unnamed" else df.columns[i] for i in range(len(df.columns))]
failer=df[df['Rang\nRank']=='RET']['Nat. / Voile\nNat. / Sail'].map(lambda x: x.split(' ')[1]).values
df=df[df.index>0].dropna(axis=1, how="all").dropna(axis=0, how="any")
return df,failer
def create_df_all_ranking_files():
df, failer=get_df_from_excel(sorted(os.listdir('downloaded_rank_files'))[::-1][0])
df["date"]=sorted(os.listdir('downloaded_rank_files'))[::-1][0][:-5]
failers=[int(f) for f in failer]
for ranking_file in sorted(os.listdir('downloaded_rank_files'))[::-1][1:-1]:
df_temp, failer=get_df_from_excel(ranking_file)
failers+=[int(f) for f in failer if int(f) not in failers]
df_temp['date']=ranking_file[:-5]
df=df.append(df_temp)
return clean_df(df, failers)
def clean_df(df, failers):
df.columns=df.columns.map(lambda x: x if x.find('\n')==-1 else x[:x.find('\n')])
df['position']=tuple(map(convert_longitude_latitude_to_float,df['Latitude'], df['Longitude']))
df['Latitude']=df['position'].map(lambda x: x[0])
df['Longitude']=df['position'].map(lambda x: x[1])
df, part_added=rename_cols(df)
df=val_to_float(df, part_added)
df["Nat"]=df['Nat. / Voile'].map(lambda x: x[x.find(' ')-3:x.find(' ')])
df["Voile"]=df['Nat. / Voile'].map(lambda x: int(x[1+x.find(' '):]))
df["Skipper"]=df["Skipper / Bateau"].map(lambda x: x[:x.find('\n')])
df["Bateau"]=df["Skipper / Bateau"].map(lambda x: x[x.find('\n'):])
df=df.drop(columns=['Nat. / Voile',"Skipper / Bateau"])
df["Rang"]=df["Rang"].map(int)
df=df[df['Voile'].map(lambda x: True if x not in failers else False)]
df=df.set_index(df["Voile"].map(str)+df["date"])
return df
def rename_cols(df):
cols=df.columns.values
part_added=[]
for col in range(len(cols)):
if "Depuis" in cols[col] :
cols[col]="Cap "+cols[col]
for col in range(len(cols)):
if cols[col]=='Vitesse':
cols[col]=cols[col]+" "+cols[col-1][4:]
part_added.append(cols[col-1][4:])
elif cols[col]=='VMG':
cols[col]=cols[col]+" "+cols[col-2][4:]
elif cols[col]=='Distance':
cols[col]=cols[col]+" "+cols[col-3][4:]
df.columns=cols
return df,part_added
def val_to_float(df, part_added):
#make the Vitesse, VMG, and distance as float
for part in part_added:
for col_type in ['Vitesse','VMG', 'Distance']:
df[col_type+' '+part]=df[col_type+' '+part].map(lambda x: float(x[: x.find(' ')]))
for col_type in ['DTF','DTL']:
df[col_type]=df[col_type].map(lambda x: float(x[: x.find(' ')]))
return df
def convert_longitude_latitude_to_float(lattitude,longitude):
#Convert longitude and lattitude in a tupple of float
return convert_latitude(lattitude),convert_longitude(longitude)
def check_if_nan_in_list(my_list):
#function that return true if there is a Nan in a list, False if not
#Not use in scrapping but later in the project
return sum([1 for item in my_list if np.isnan(item)])>0
def convert_latitude(latitude):
#formula to convert latitude in float
d, m, s = map(float, (latitude[:2], latitude[3:5], latitude[6:8]))
return (d + m / 60. + s / 3600.) * (1 if ('N' in latitude) else -1)
def convert_longitude(longitude):
#formula to convert longitude in float
d, m, s = map(float, (longitude[:2], longitude[3:5], longitude[6:8]))
return (d + m / 60. + s / 3600.) * (1 if ('E' in longitude) else -1)
def add_boat_infos(df, boat_dict):
#Function that create a dataframe by merging the dataframe created with ranking files and the dictionnary from the glossaire page
df_complete=copy.deepcopy(df)
for key in boat_dict.get(list(boat_dict.keys())[np.argmin(list(map(lambda x: len(boat_dict.get(x)), boat_dict.keys())))]).keys():
df_complete[key]=df_complete['Voile'].map(lambda x: boat_dict.get(x).get(key))
df_complete['Année lancement']=df_complete['Date de lancement'].map(lambda x: int(x[-4:]))
df_complete['Longueur']=df_complete['Longueur'].map(lambda x: float(x[: x.find(' ')].replace(',','.')))
df_complete['Largeur']=df_complete['Largeur'].map(lambda x: float(x[: x.find(' ')].replace(',','.')))
df_complete["Tirant d'eau"]=df_complete["Tirant d'eau"].map(lambda x: float(x[: x.find(' ')].replace(',','.')))
df_complete['Déplacement (poids)']=df_complete['Déplacement (poids)'].map(lambda x: float(x[: x.find(' ')].replace(',','.')) if x!='NC' and x!='nc' else np.nan)
df_complete['Hauteur mât']=df_complete['Hauteur mât'].map(lambda x: float(x[: x.find(' ')].replace(',','.')))
df_complete['Surface de voiles au près']=df_complete['Surface de voiles au près'].map(lambda x: float(x[: x.find(' ')].replace(',','.')))
df_complete['Surface de voiles au portant']=df_complete['Surface de voiles au portant'].map(lambda x: float(x[: x.find(' ')].replace(',','.')))
df_complete['Foil']=df_complete['Nombre de dérives'].map(lambda x: True if x in ['foils', 'foiler'] else False)
return df_complete
Récupération des données et mise en forme dans le dataframe¶
if download_last_datas:
download_files()
#Creation of a dataframe using the ranking file
df=create_df_all_ranking_files()
#Get the specific informations about the boats on the glossaire page
boat_dict=request_boats_info(glossaire_url)
#Merge all those information in a dataframe
df_complete=add_boat_infos(df, boat_dict)
Partie III - Utilisation des données au travers de modèles¶
A- Utilisation des données du bateau pour estimer et sa VMG¶
#We will try a regression beetween the informations we have about the boat, in our dataframe we currently have this columns
print(df_complete.columns)
# In order to don't have overlapping beetween speed of different classement we will
#try to explain the ones "Depuis le dernier classement"
#We can use the informations Longueur, Largeur, Tirant d'eau, Année de lancement, Surface de voile au près et au portant
#Nombre de dérives, Hauteur Mât, and deplacement poids
#However we can see that Nombre de dérives doesn't only includes int, but also "foils" or "foiler"
print("distinct values in 'Nombre de dérives' "+str(df_complete['Nombre de dérives'].unique()))
#A solution to this is to create a dummies variables
#But as foils and foiler refer to the same things, we will use the created column Foil that contain true of false if the boat have
#foils or not
print("distinct values in 'Foil' "+str(df_complete['Foil'].unique()))
#We can also see that some informations are missing in Déplacement poids:
print("distinct values in 'Déplacement (poids)' "+str(df_complete['Déplacement (poids)'].unique()))
print("There is nan in Déplacement (poids): "+ str(check_if_nan_in_list(list(df_complete['Déplacement (poids)'].unique()))))
#A solution to this is to change the nan by the mean values of this columns
#Let's first create the explained variables of our two models
Y_VMG=df_complete['VMG Depuis le dernier classement'].values.reshape(-1,1)
Y_Vitesse=df_complete['Vitesse Depuis le dernier classement'].values.reshape(-1,1)
#Let's check if both don't contain NaN values
print("There is nan in Y_VMG: "+str(check_if_nan_in_list(Y_VMG)))
print("There is nan in Y_Vitesse: "+str(check_if_nan_in_list(Y_Vitesse)))
#We can see that we don't have nan, we can go to the explanatory variables
#Transform the foil information in a dummy variable
foils=df_complete['Foil'].map(lambda x: 1 if x else 0).values
#get the weights of boats and change nan with the mean
mean_weight=np.mean([item for item in df_complete['Déplacement (poids)'].unique() if np.isnan(item)==False])
weights=df_complete['Déplacement (poids)'].map(lambda x: mean_weight if np.isnan(x) else x).values
X=np.zeros((len(weights),9 ))
X[:,-1]=foils
X[:,0]=weights
X[:,1:-1]=df_complete[['Longueur', 'Largeur', "Tirant d'eau", 'Hauteur mât',
'Surface de voiles au près', 'Surface de voiles au portant',
'Année lancement']]
X_colnames=['Déplacement (poids)','Longueur', 'Largeur', "Tirant d'eau", 'Hauteur mât',
'Surface de voiles au près', 'Surface de voiles au portant',
'Année lancement','Foil']
clf=linear_model.LinearRegression(fit_intercept=True, normalize=True)
clf.fit(X,Y_VMG)
print("Les coefficients de la regressions de la VMG par rapport à nos variables sont:")
print([X_colnames[i]+": "+str(clf.coef_[0,i])+", " for i in range(len(X_colnames))])
res=Y_VMG-clf.predict(X)
plt.plot(res)
plt.xlabel("observation number")
plt.ylabel("Residus")
plt.title("Resultats avec une regression linéaire classique")
plt.show()
lasso_model = LassoCV(cv=5, random_state=0,n_alphas=100, normalize=True).fit(X, Y_VMG.reshape(-1))
variable_selected_lasso=[X_colnames[var]+": "+str(lasso_model.coef_[var]) for var in range(len(lasso_model.coef_)) if np.isclose(lasso_model.coef_[var],0)==False]
residus_test_lasso=lasso_model.predict(X)-Y_VMG.reshape(-1)
print("Les coefficients de la regressions de la VMG en utilisant un modèle Lasso par rapport à nos variables sont:")
print(variable_selected_lasso)
plt.plot(residus_test_lasso)
plt.xlabel("observation number")
plt.ylabel("Residus")
plt.title("Resultats avec une regression lasso")
plt.show()
En regardant les résidus on peut voir que ceux-ci sont assez importants, mais surtout que ceux-ci semblent heteroscedastiques. En effet, de par la création de notre dataframe les résidus se trouvent être par date. Hors la vitesse d'un bateau ne dépend pas uniquement de ses caracteristiques mais aussi des conditions météorologiques ! Cette première regression nous permet donc d'estimer une VMG moyenne du bateau, et l'utilisation d'un modèle Lasso de trouver les variables significatives pour determiner cela. Cela ne signifie pas que les autres variables n'entrent pas en compte dans la réalité, toutefois les caracteristiques étant corrélées entre elles, prendre toutes celles-ci n'a pas d'interet.
Dans un premier temps nous pouvons donc estimer la VMG moyenne que devrait avoir chacun des bateaux en fonction de ses caracteristiques:
dt=df_complete['date'].max()
df_stats_boats=df_complete[df_complete['date']==dt]
foils=df_stats_boats['Foil'].map(lambda x: 1 if x else 0).values
#get the weights of boats and change nan with the mean
mean_weight=np.mean([item for item in df_stats_boats['Déplacement (poids)'].unique() if np.isnan(item)==False])
weights=df_stats_boats['Déplacement (poids)'].map(lambda x: mean_weight if np.isnan(x) else x).values
X=np.zeros((len(weights),9 ))
X[:,-1]=foils
X[:,0]=weights
X[:,1:-1]=df_stats_boats[['Longueur', 'Largeur', "Tirant d'eau", 'Hauteur mât',
'Surface de voiles au près', 'Surface de voiles au portant',
'Année lancement']]
VMG_boats=lasso_model.predict(X)
df_Expected_VMG=pd.DataFrame({'Voile':df_stats_boats['Voile'].values,'Rang': df_stats_boats['Rang'].values, 'VMG Expected': VMG_boats })
boats_id=df_Expected_VMG['Voile']
Nous pouvons désormais voir quel bateau réussi à battre sa VMG expected (bon choix du skipper), et ceux qui font moins bien (choix moins judicieux du skipper)
realized_mean_VMG=[df_complete[df_complete['Voile']==boats_id[id]]['VMG Depuis le dernier classement'].mean() for id in range(len(boats_id))]
df_Expected_realized_VMG=pd.DataFrame({'Voile':df_stats_boats['Voile'].values,'Rang': df_stats_boats['Rang'].values, 'VMG Expected': VMG_boats, 'VMG realized': realized_mean_VMG })
df_Expected_realized_VMG['Difference Realized Expected']=df_Expected_realized_VMG['VMG realized']-df_Expected_realized_VMG['VMG Expected']
print(df_Expected_realized_VMG.sort_values(by='Difference Realized Expected', ascending=False))
En enlevant l'effet du bateau sur la VMG on peut voir que le classement ne serait alors plus le même. Cette transformation a l'avantage de permettre une comparaison plus équitable entre les skippers. Cependant nous n'en sommes encore qu'au début du Vendée Globe, hors une VMG plus faible n'est pas necessairement signe d'un moins bon skipper à l'heure actuelle, le plus important étant la VMG moyenne jusqu'à l'arrivée. En effet, il est possible et probable que certains skippers aient fait le choix d'une VMG plus faible, par exemple en s'éloignant du tracé le plus court, pour aller prendre des vents plus puissants par la suite.
Observons désormais la relation entre la difference entre réalisée et attendue entre les deux parties de la course
dt_sorted=np.sort(df_complete['date'].unique())
dt_first_part=dt_sorted[:int(len(dt_sorted)/2)]
dt_last_part=dt_sorted[int(len(dt_sorted)/2):]
df_complete_first_part=df_complete[df_complete['date'].map(lambda x: True if x in dt_first_part else False) ]
df_complete_last_part=df_complete[df_complete['date'].map(lambda x: True if x in dt_last_part else False) ]
realized_mean_VMG_first_part=[df_complete_first_part[df_complete_first_part['Voile']==boats_id[id]]['VMG Depuis le dernier classement'].mean() for id in range(len(boats_id))]
realized_mean_VMG_last_part=[df_complete_last_part[df_complete_last_part['Voile']==boats_id[id]]['VMG Depuis le dernier classement'].mean() for id in range(len(boats_id))]
#df_Expected_realized_VMG_last_part=pd.DataFrame({'Voile':df_stats_boats['Voile'].values,'Rang': df_stats_boats['Rang'].values, 'VMG Expected': VMG_boats last part, 'VMG realized first part': realized_mean_VMG_first_part,'VMG realized last part': realized_mean_VMG_last_part })
df_Expected_realized_VMG_by_part=pd.DataFrame({'Voile':df_stats_boats['Voile'].values,'Rang': df_stats_boats['Rang'].values, 'VMG Expected': VMG_boats, 'VMG realized first part': realized_mean_VMG_first_part, 'VMG realized last part': realized_mean_VMG_last_part })
df_Expected_realized_VMG_by_part['Difference Realized Expected first part']=df_Expected_realized_VMG_by_part['VMG realized first part']-df_Expected_realized_VMG_by_part['VMG Expected']
df_Expected_realized_VMG_by_part['Difference Realized Expected last part']=df_Expected_realized_VMG_by_part['VMG realized last part']-df_Expected_realized_VMG_by_part['VMG Expected']
df_Expected_realized_VMG_by_part['Difference Realized Expected variation']=df_Expected_realized_VMG_by_part['Difference Realized Expected last part']-df_Expected_realized_VMG_by_part['Difference Realized Expected first part']
Y=df_Expected_realized_VMG_by_part['Difference Realized Expected variation'].values.reshape(-1,1)
X=df_Expected_realized_VMG_by_part['Difference Realized Expected first part'].values.reshape(-1,1)
clf=linear_model.LinearRegression(fit_intercept=True)
clf.fit(X,Y)
df_Expected_realized_VMG_by_part
print("intercept is: "+str(clf.intercept_)+" and coefficient is: " +str(clf.coef_))
On peut voir sur la seconde partie de la course que l'ensemble des bateaux affichent une VMG supérieure (du fait des conditions météorologiques dans la zone où ceux-ci se trouvent). Cet effet se retrouve dans l'intercept de notre regression. Ce qu'il est interressant de voir aussi est l'effet de rattrapage, en effet, le coefficient de notre régression est négatif, cela veut dire que les bateaux ayant été le plus lentement par rapport à l'attendue sur la premiere partie sont ceux qui ont été les plus rapides sur la seconde partie de la course. On peut donc en déduire que l'hypothèse faite precedement se vérifie.
Cependant cette analyse souffre d'un biais important. Utiliser en premier lieu la VMG qui dépend du cap, et donc des choix du skipper ne donne pas un indicateur fiable de la performance des bateaux, car elle inclut fortement les choix du skipper. Afin de corriger ce biais j'ai choisi d'essayer d'estimer par bateaux la vitesse en fonction du vent et de son sens. Cela permet une comparaison plus fine des performances des bateaux en dehors des performances humaines. Ici la difficulté vient bien plus de la récupération de la donnée, et du raccrochement à nos données initiales. Le code ci-dessous se charge de recupérer les données et de les ajouter à notre dataframe. Les données proviennent du NOAA (National Oceanic and Atmospheric Administration), de deux fichiers en format .nc dont l'un contient la vitesse du vent méridionnale et le second du vent zonal. Ces données sont observées 4x par jour sur une grille de longitude et de lattitude. De ces deux valeurs en m.s-1 est recalculé un angle du vent dans le même référentiel que celui utilisé pour le cap du bateau, ainsi que sa force. 3 colonnes sont alors ajoutées contenant l'angle du vent, l'angle du vent par rapport au bateau et sa vitesse, et ce à chaque date en fonction des données les plus proches géographiquement et temporellement.
B- Utilisation des données de vents pour estimer la VMG¶
def download_wind_data_files():
url_current_year_Vvalues="ftp://ftp2.psl.noaa.gov/Datasets/ncep.reanalysis/surface_gauss/vwnd.10m.gauss.2020.nc"
url_last_year_Vvalues="ftp://ftp2.psl.noaa.gov/Datasets/ncep.reanalysis/surface_gauss/vwnd.10m.gauss.2019.nc"
url_current_year_Uvalues="ftp://ftp2.psl.noaa.gov/Datasets/ncep.reanalysis/surface_gauss/uwnd.10m.gauss.2020.nc"
url_last_year_Uvalues="ftp://ftp2.psl.noaa.gov/Datasets/ncep.reanalysis/surface_gauss/uwnd.10m.gauss.2019.nc"
delete_if_exist_and_recreate_directory(os.getcwd()+"\\wind_datas\\")
urllib.request.urlretrieve(url_current_year_Vvalues, os.getcwd()+"\\wind_datas\\current_year_Vvalues.nc")
urllib.request.urlretrieve(url_last_year_Vvalues, os.getcwd()+"\\wind_datas\\last_year_Vvalues.nc")
urllib.request.urlretrieve(url_current_year_Uvalues, os.getcwd()+"\\wind_datas\\current_year_Uvalues.nc")
urllib.request.urlretrieve(url_last_year_Uvalues, os.getcwd()+"\\wind_datas\\last_year_Uvalues.nc")
def get_lon_lat_time_from_nc(nc_file):
#function returning the list of available longitude and latitude in the file
variables_in_files=nc_file.variables
return variables_in_files['lon'][:].data,variables_in_files['lat'][:].data,variables_in_files["time"][:].data
def get_values_wind(Unc_file, Vnc_file, lat_in_file, lon_in_file, time_in_file, real_lat, real_lon, real_time ):
#Using a real longitude, latitude and time wanted, return the values the closest in the space and time from the nc file
lat_index=get_closest_index(lat_in_file, real_lat)
if real_lon>0:
real_lon=real_lon
else:
real_lon=360+real_lon
lon_index=get_closest_index(lon_in_file, real_lon)
time_index=get_closest_index(time_in_file, real_time)
return calculate_wind_speed_angle(Unc_file.variables['uwnd'][time_index, lat_index,lon_index].data, Vnc_file.variables['vwnd'][time_index,lat_index,lon_index].data)
def calculate_wind_speed_angle(Uwind, Vwind):
#Use zonal and meridional wind to calculate a wind speed and angle
wind_speed=np.sqrt(Uwind**2+Vwind**2)
if wind_speed!=0:
wind_angle=put_360(math.degrees(np.arccos(-Vwind/wind_speed)))
else:
wind_angle=0
return wind_speed, wind_angle
def put_360(degrees_180):
#convertion function
if degrees_180<0:
degrees_180=360-degrees_180
return degrees_180
def get_closest_index(available_values, searched_value):
#return the index of the closest value in a list of a searched value
return np.argmin(abs(available_values-searched_value))
def transform_date(date_df):
y=int(date_df[:4])
m=int(date_df[4:6])
d=int(date_df[6:8])
h=int(date_df[9:11])
return (datetime(year=y, month=m, day=d, hour=h)-datetime(2020,1,1)).total_seconds()/3600
if download_last_datas:
download_wind_data_files()
cy_Vwind = Dataset(os.getcwd()+"\\wind_datas\\current_year_Vvalues.nc", "r", format="NETCDF4")
cy_Uwind = Dataset(os.getcwd()+"\\wind_datas\\current_year_Uvalues.nc", "r", format="NETCDF4")
ly_Vwind = Dataset(os.getcwd()+"\\wind_datas\\last_year_Vvalues.nc", "r", format="NETCDF4")
ly_Uwind = Dataset(os.getcwd()+"\\wind_datas\\last_year_Uvalues.nc", "r", format="NETCDF4")
lon_in_file,lat_in_file, time_in_file_cy=get_lon_lat_time_from_nc(cy_Vwind)
time_in_file_cy-=time_in_file_cy[0]
_, _, time_in_file_ly=get_lon_lat_time_from_nc(ly_Uwind)
time_in_file_ly-=time_in_file_ly[0]
df_complete['time_delta_start_year']=df_complete['date'].map(transform_date)
real_lats=df_complete['Latitude'].values
real_lons=df_complete['Longitude'].values
real_times=df_complete['time_delta_start_year'].values
wind_speed=[]
wind_angle=[]
for iter in range(len(real_times)):
vspeed, vangle=get_values_wind(cy_Uwind,cy_Vwind, lat_in_file, lon_in_file, time_in_file_cy, real_lats[iter],real_lons[iter],real_times[iter])
wind_speed.append(vspeed)
wind_angle.append(vangle)
df_complete['Wind speed']=wind_speed
df_complete['wind angle']=wind_angle
df_complete['wind angle with boat']=df_complete['wind angle']-df_complete['Cap Depuis le dernier classement'].map(lambda x: int(x[:-1]))
df_complete['wind angle with boat']=df_complete['wind angle with boat'].map(lambda x: abs(x) if abs(x)<180 else abs(x)-90).map(lambda x: x if x<180 else abs(x-180)).map(lambda x: abs(x-180))
print(df_complete[df_complete['Rang']==1][df_complete['date']=="20201116_170000"]['Wind speed'])
print(df_complete[df_complete['Rang']==1][df_complete['date']=="20201116_170000"]['wind angle'])
df_complete
La première approche, intuitive serait de regresser uniquement la vitesse et l'angle du vent sur le bateau. Cependant les résultats sont évidement assez mauvais, car la force appliquée par le vent n'est pas linéaire de ces deux variables. Un dataframe temporaire est ici créé avec des variables tel que le sinus et le cosinus de l'angle du vent. Enfin le bateau n'avancant pas sans vent, contrairement à la plupart des modèles il faut ici bien enlever l'intercept qui n'aurait alors pas de sens. Enfin les données sont réduites en séparant le dataframe par bateau (100 points à l'heure actuelle), et les skippers ne se mettant pas face au vent certaines configurations ne sont pas représentées ou dans des zones de vent nul. Pour contrer cela, j'ai agrandi les données d'étude en y ajoutant des valeurs dans ces cas-ci, afin de forcer le modèle à respecter ces contraintes. Enfin pour utiliser les interactions entre les variables et les variables au carré j'ai fait appel à PolynomialFeatures de sklearn.
from sklearn.preprocessing import PolynomialFeatures
boats_numbers=df_complete['Voile'].unique()
df_for_reg=copy.deepcopy(df_complete)
dt_part=dt_sorted[3:-3]
df_for_reg=df_for_reg[df_for_reg['date'].map(lambda x: True if x in dt_part else False)]
df_for_reg['div']=df_for_reg['wind angle with boat'].map(lambda x: 1/np.log(x+2))
df_for_reg['sin']=df_for_reg['wind angle with boat'].map(lambda x: np.sin(math.radians(x)))
df_for_reg['cos']=df_for_reg['wind angle with boat'].map(lambda x: np.cos(math.radians(x)))
clf=linear_model.LinearRegression(fit_intercept=False)
R2=[]
boats_model_lasso={}
boats_model_clf={}
R22=[]
test=[]
poly = PolynomialFeatures(degree=2,include_bias=False, interaction_only=False)
for boat in boats_numbers:
Y=df_for_reg[df_for_reg['Voile']==boat][df_for_reg['Vitesse Depuis le dernier classement']>0.01]['Vitesse Depuis le dernier classement']
X=df_for_reg[df_for_reg['Voile']==boat][df_for_reg['Vitesse Depuis le dernier classement']>0.01][['wind angle with boat', 'Wind speed', 'div','cos','sin']]
Xbis=np.zeros((int(np.shape(X)[0]*1.3),np.shape(X)[1]))
Xbis[:np.shape(X)[0],:]=X
Xbis[np.shape(X)[0]:,1]=0
Xbis[np.shape(X)[0]:,0]=np.random.randint(0,180,np.shape(Xbis)[0]-np.shape(X)[0])
Xbis[np.shape(X)[0]:,2]=1/np.log(np.random.randint(0,180,np.shape(Xbis)[0]-np.shape(X)[0])+2)
Xbis[np.shape(X)[0]:,3]=np.sin(np.radians(np.random.randint(0,180,np.shape(Xbis)[0]-np.shape(X)[0])))
Xbis[np.shape(X)[0]:,4]=np.cos(np.radians(np.random.randint(0,180,np.shape(Xbis)[0]-np.shape(X)[0])))
Ybis=np.zeros((len(Xbis)))
Ybis[:np.shape(X)[0]]=Y
Ybis[np.shape(X)[0]:]=0
#Xbis=X.values
#Ybis=Y.values
Xbis=poly.fit_transform(Xbis)
lasso_model=LassoCV(fit_intercept=False, normalize=True).fit(Xbis, Ybis)
#lasso_model = linear_model.Lasso(max_iter=100000, selection='random',fit_intercept=False, normalize=False).fit(Xbis, Ybis)
clf=linear_model.LinearRegression(fit_intercept=False).fit(Xbis,Ybis)
pred2=lasso_model.predict(Xbis)
pred=clf.predict(Xbis)
#res=Y-pred
#res2=Y-pred2
plt.plot(range(len(pred)),pred, label='predit')
plt.plot(range(len(pred)),pred2, label='predit Lasso')
plt.plot(range(len(pred)),Ybis, label='True')
R2.append(metrics.r2_score(Ybis, pred))
R22.append(metrics.r2_score(Ybis, pred2))
plt.legend()
plt.title("Valeurs réelle vs estimé pour le bateau numéro "+str(boat))
plt.show()
boats_model_clf.update({boat: clf})
boats_model_lasso.update({boat: lasso_model})
test.append((np.min(X.values[:,1]),np.max(X.values[:,1]), np.mean(X.values[:,1])))
print("R2 moyen obtenu: "+str(np.mean(R2)))
print("R2 moyen obtenu: "+str(np.mean(R22)))
On peut voir que les résultats de la regression, sans être totalement parfaits semblent satisfaisants. Il était dur d'esperer beaucoup plus étant donné d'une part l'approximation des données au niveau géographique et temporel, mais aussi d'autres données météorologiques pouvant ralentir la navigation comme la houle, les précipitations, ou les changements de cap opérés par le skipper. Cependant ces résultats me semblent suffisants pour pouvoir comparer les bateaux entre eux. Une première question qui nous était proposée était d'estimer l'effet des foils sur les bateaux. Je ne l'avais pas fait pour la VMG, car cela n'aurait pas eu beaucoup de sens, mais ici cela devient bien plus interressant. Je vais donc entamer la comparaison entre les bateaux pour des conditions de vents différentes. On voit ici qu'il existe tout de même une erreur importante notamement sur les zones de vent nul, avec des erreurs assez importantes.
Pour ce faire il est necessaire de prendre en compte les autres effets sur le bateaux en plus des foils. Nous allons donc regresser la vitesse prédite pour chaque bateau et dans des conditions de vents differentes sur leurs caractéristiques.
def transform_data_for_model(ws, wa):
watransform=np.log(wa+2)
return np.array([wa, ws,1/watransform, np.sin(math.radians(wa)), np.cos(math.radians(ws))]).reshape(1,-1)
df_infos_per_boat_only=df_for_reg[df_for_reg['date'].map(lambda x: True if x in dt_part[-1] else False)]
df_infos_per_boat_only=df_infos_per_boat_only[['Voile','Longueur', 'Largeur', "Tirant d'eau", 'Hauteur mât',
'Surface de voiles au près', 'Surface de voiles au portant',
'Année lancement','Foil']]
nb_boat_with_foil=df_infos_per_boat_only['Foil'].sum()
nb_boat_without_foil=len(df_infos_per_boat_only)-df_infos_per_boat_only['Foil'].sum()
wind_speed_step=0.25
wind_speed_max=18
wind_speed_min=0
wind_angles=[0,45,90,135,180]
speed_boat_with_foil=np.zeros((int((wind_speed_max-wind_speed_min)/wind_speed_step), len(wind_angles)))
speed_boat_without_foil=np.zeros((int((wind_speed_max-wind_speed_min)/wind_speed_step), len(wind_angles)))
Y=np.zeros((int((wind_speed_max-wind_speed_min)/wind_speed_step*len(boats_numbers))))
X=np.zeros((int((wind_speed_max-wind_speed_min)/wind_speed_step*len(boats_numbers)), 8))
windangle=90
col_names=['Longueur', 'Largeur', "Tirant d'eau", 'Hauteur mât',
'Surface de voiles au près', 'Surface de voiles au portant',
'Année lancement','Foil']
col=0
for windangle in wind_angles:
count=0
for boat in boats_numbers:
clf=boats_model_lasso.get(boat)
line=0
for windspeed in np.arange(wind_speed_min, wind_speed_max, wind_speed_step):
Y[count]=clf.predict(poly.fit_transform([transform_data_for_model(windspeed, windangle)[0]]))
X[count,:]=df_infos_per_boat_only[df_infos_per_boat_only['Voile']==boat].values[0,1:]
if X[count,-1]:
speed_boat_with_foil[line, col]+=Y[count]/nb_boat_with_foil
else:
speed_boat_without_foil[line, col]+=Y[count]/nb_boat_without_foil
line+=1
count+=1
lasso_model = LassoCV(cv=5, random_state=0,n_alphas=100, fit_intercept=False, normalize=True).fit(X, Y)
print("les variables du modeles lasso pour un angle de "+str(windangle)+" degré sont :")
print(lasso_model.coef_)
print("les variables selectionnés sont ainsi:")
print([col_names[var] for var in range(len(lasso_model.coef_)) if np.isclose(lasso_model.coef_[var],0)==False])
col+=1
for col in range(len(wind_angles)):
plt.plot(np.arange(wind_speed_min, wind_speed_max, wind_speed_step),speed_boat_with_foil[:,col], label='bateaux avec foil')
plt.plot(np.arange(wind_speed_min, wind_speed_max, wind_speed_step),speed_boat_without_foil[:,col], label='bateaux sans foil')
plt.legend()
plt.xlabel("Vitesse du vent (m/s)")
plt.ylabel("Vitesse du navire (noeuds)")
plt.title("resultat avec un angle d'incidence du vent de "+str(wind_angles[col])+" degrés")
plt.show()
Bien que le modèle soit loin de présenter une précision absolue dans le sens où de nombreuses conditions sont ignorées, on peut voir que les résultats montrent tout de même une certaine cohérence sur certains points. Les courbes partent toutes d'environ 0, en vent nul le bateau n'avance donc presque pas. Quand celui-ci augmente sa vitesse augmente et ce en fonction de l'angle du vent. Face au vent de même le bateau n'avance presque pas. Enfin il est interressant que bien qu'un lasso ne détecte pas l'effet des foils directement mais l'année de lancement du bateau comme explicatif de la vitesse, l'analyse graphique des courbes de vitesses prédites permet de voir que les bateaux avec foils sont légèrement avantagés sur les autres. On peut d'ailleurs voir que les bateaux possédant des foils sont principalement les bateaux récents. L'information année de lancement contient donc en partie l'information foil, mais aussi d'autres informations expliquant sa séléction à contrario de la variable foil.
Lien entre Vitesse et VMG¶
Y=df_complete['VMG Depuis le dernier classement'].values.reshape(-1,1)
X=df_complete["Vitesse Depuis le dernier classement"].values.reshape(-1,1)
clf=linear_model.LinearRegression(fit_intercept=False).fit(X,Y)
pred=clf.predict(X)
plt.plot(Y, label='True')
plt.plot(pred, label='predit')
plt.legend()
plt.show()
print("Intercept: "+str(clf.intercept_))
print("coefficient "+str(clf.coef_))
print("R2 de la regression: "+str(metrics.r2_score(pred,Y)))
Le lien entre vitesse et VMG est assez évident mais il apparait tout de même utile de le verifier, on voit ici que le modèle fonctionne très bien avec un R2 important. Ces deux variables auraient d'ailleurs pu être moins liées dans le cas où les routes prises par les skippers auraient fortement différé, auquel cas une vitesse plus rapide aurait pu entrainer une VMG plus petite car fesant un détour important. Ce n'est ici pas le cas, cela nous permet donc de voir que la vitesse des bateaux est ici un élément déterminant dans la VMG. La VMG impliquant méchaniquement le classement, l'observation de la vitesse théorique de chaque bateau a donc bien eu un interêt certain pour comprendre le classement.
Ayant désormais le lien entre la VMG et la vitesse, et ayant un modèle, bien qu'imparfait, de calcul théorique de la vitesse pour chaque bateau il est possible de recréer la VMG attendue pour chaque bateau dans des conditions similaires à celles rencontrées sur la course. Pour ce faire j'ai récupéré les données de vents de tous les bateaux durant la course. La fonction get_expected_vmg_boat calcule alors la VMG du bateau pour ces données en fonction de la vitesse prédite de celui-ci, et en fait la moyenne pour l'ensemble des données de vent.
wind_speed_datas=df_complete['Wind speed'].values
wind_angle_datas=df_complete['wind angle with boat'].values
def get_expected_vmg_boat(clf_vmg_vitesse, boat_speed_model, winds_speed_datas, wind_angle_datas):
vmg_e=[]
for iter in range(len(wind_speed_datas)):
vmg_e.append(clf_vmg_vitesse.predict(boat_speed_model.predict(poly.fit_transform([transform_data_for_model(winds_speed_datas[iter], wind_angle_datas[iter])[0]])).reshape(-1,1))[0][0])
return np.mean(vmg_e)
Expected_VMG=[]
for boat in boats_numbers:
#Expected_VMG.append(clf.predict(boats_model.get(boat).predict(poly.fit_transform([transform_data_for_model(mean_wind, mean_wind_angle)[0]])).reshape(-1,1))[0][0])
Expected_VMG.append(get_expected_vmg_boat(clf, boats_model_lasso.get(boat), wind_speed_datas,wind_angle_datas))
df_Expected_VMG=pd.DataFrame({'Voile':df_stats_boats['Voile'].values,'Rang': df_stats_boats['Rang'].values, 'VMG Expected': VMG_boats, 'New VMG expected': Expected_VMG })
realized_mean_VMG=[df_complete[df_complete['Voile']==boats_id[id]]['VMG Depuis le dernier classement'].mean() for id in range(len(boats_id))]
Foil=[df_complete[df_complete['Voile']==boats_id[id]]['Foil'].values[0] for id in range(len(boats_id))]
Annee=[df_complete[df_complete['Voile']==boats_id[id]]['Année lancement'].values[0] for id in range(len(boats_id))]
df_Expected_realized_VMG=pd.DataFrame({'Voile':df_stats_boats['Voile'].values,'Rang': df_stats_boats['Rang'].values, 'VMG Expected': VMG_boats, 'New VMG expected': Expected_VMG, 'VMG realized': realized_mean_VMG, 'Foil': Foil, 'Année lancement':Annee })
df_Expected_realized_VMG['Difference Realized Expected']=df_Expected_realized_VMG['VMG realized']-df_Expected_realized_VMG['VMG Expected']
df_Expected_realized_VMG['Difference Realized New Expected']=df_Expected_realized_VMG['VMG realized']-df_Expected_realized_VMG['New VMG expected']
df_Expected_realized_VMG['Difference Expected New Expected']=df_Expected_realized_VMG['VMG Expected']-df_Expected_realized_VMG['New VMG expected']
df_Expected_realized_VMG
En recréant le premier tableau on voit ici que la différence entre le réalisé présentent un biais haussier, dû en parti au fait que le modèle n'inclut pas d'intercept ici contrairement au premier. Le modele n'est donc pas parfait, durant mes recherches d'ailleurs en utilisant d'autres combinaisons de variables, et en utilisant des intercepts les résultats ici semblaient plus interressant dans le sens où ils étaient en absolu plus proche des valeurs réalisé et moins biaiser. Cependant ces modèle présentaient alors moins de sens car présentant des vitesses pouvant être positive en vent nul, et pouvant passer négative en cas de vent trop important.
print("sur l'ensemble des bateaux")
print("moyenne realisé: "+str(df_Expected_realized_VMG['VMG realized'].mean()))
print("moyenne expected: "+str(df_Expected_realized_VMG['VMG Expected'].mean()))
print("moyenne new expected: "+str(df_Expected_realized_VMG['New VMG expected'].mean()))
En regardant sur l'ensemble des bateaux, evidemment la premiere facon de calculer la VMG est par construction égale à celle réalisé (regression avec intercept), et différentes pour celle que nous vennons de créer.
print("sur les 10 premiers")
print("moyenne realisé: "+str(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']<10]['VMG realized'].mean()))
print("moyenne expected: "+str(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']<10]['VMG Expected'].mean()))
print("moyenne new expected: "+str(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']<10]['New VMG expected'].mean()))
On voit que l'en ne regardant plus que les 10 premiers bateaux la nouvelle formule permet de calculer une VMG semblent assez similaire aux valeurs obtenu à l'aide de l'ancienne formule.
print("sur les 10 derniers (hors bateau 8)")
print("moyenne realisé: "+str(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']>21][df_Expected_realized_VMG['Rang']<32]['VMG realized'].mean()))
print("moyenne expected: "+str(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']>21][df_Expected_realized_VMG['Rang']<32]['VMG Expected'].mean()))
print("moyenne new expected: "+str(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']>21][df_Expected_realized_VMG['Rang']<32]['New VMG expected'].mean()))
Sur la fin de peloton toutefois les valeurs sont bien plus faible qu'avec l'ancienne formule. L'avantage du modele pour comparer les competences est ici clair, il n'y a plus de biais en fonction de la VMG réalisé appliqué par le modèle qui pénalise mécaniquement les derniers.
plt.plot(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']<32]['Difference Realized Expected'], label='ecart à expected')
plt.plot(df_Expected_realized_VMG[df_Expected_realized_VMG['Rang']<32]['Difference Realized New Expected'], label='ecart à new expected')
plt.legend()
plt.show()
En tracant le graphe des écarts avec les deux méthode la différence est frappante. En effet, là où le premier modele n'était fitter en réalité qu'en moyenne et présentait des écarts importants pour les valeurs excentré, ce nouveau modèle présente des écarts en valeurs plus constant, même si présentant un biais.
df_Expected_realized_VMG.sort_values(by='Difference Realized New Expected', ascending=False)[['Voile', 'Rang', 'VMG realized','Difference Realized New Expected', 'Foil','Année lancement']].reset_index(drop=True)
Ce tableau nous permet ainsi d'obtenir un nouveau classement, qui permettrait d'être plus juste en fonction des bateaux, et donc de comparer les compétences des skippers.
Pour conclure, ce genre de méthode permet en effet de bien rendre compte de l'impact des bateaux dans le classement de la course, ce qui n'est pas très surprennant. Il pourrait ici paraitre attrayant de corriger comme il a été fait dans ce projet les classements en fonction des performances des bateaux. Toutefois il se pose plusieurs limites à cela. D'une part la fiabilité de ces modèles, en effet dans le cas présent étant donné du manque de données (seulement 90 points par bateaux actuellement) le modèle représentant la vitesse de ceux-ci en fonction du vent reste peu fiable et très imparfait. De plus parmi les compétences des skippers le choix et/ou la création de leur navire est une composante à prendre en compte. En effet ils sont, avec leurs sponsors évidemment, les décideurs de leur navire et s'impliquent souvent dans le choix, la création et les modifications de ceux-ci.
C - Utilisation des modèles créés précédement et des données de vent pour prédire le meilleur trajet pour chaque bateau¶
Ayant désormais des modèles, bien qu'imparfaits, permettant d'estimer la vitesse des bateaux selon l'angle et la vitesse du vent, il est désormais possible d'estimer la route la plus rapide pour finir la course. Pour ce faire j'ai créé les fonctions ci-dessous. Le fonctionnement est assez simple, ayant mis des points balises sur la carte, une fonction va chercher le chemin le plus rapide d'un point à l'autre. Pour ce faire deux méthodes sont utilisées ici. La premiere est la fonction from_a_to_b qui va chercher à estimer le chemin le plus rapide en regardant toutes les 6h quel cap permet de se rapprocher le plus du point objectif en fonction des conditions de vent, en essayant un nombre de déviations défini. La seconde from_a_to_b_with_random_variation, va découper le parcours entre les balises en sous parcours en ajoutant des points incluant des déviations plus ou moins importantes. Cette méthode est appliquée de facon récursive et fait au final appel à la méthode from_a_to_b puis compare les longueurs totale des chemins pour finalement decider quel chemin est le plus rapide. Pour faire tout cela il faut évidement disposer de données de vent, hors les valeurs futures ne sont bien évidemment pas disponibles. Pour gérer cela, bien que de façon imparfaite nous utiliserons les données de vent à la même date et heure mais à l'année passée.
global cap_bonne_esperance_pos
global Leuuwin_pos
global canary
global portugal
global horn
global goal
cap_bonne_esperance_pos=(-39,17.5)
Leuuwin_pos=(-45,128)
horn=(-58,-65)
canary=(20,-26)
portugal=(44,-10)
goal=(46.26,-2.5)
def find_best_way_to_goal(cy_Uwind,cy_Vwind,ly_Uwind,ly_Vwind,start_date,boats_model, boatnumber,
nb_random_mid_point, start_post, next_point=0,nb_dev=1):
#main function to search the path from the boat to the goal
boat_speed_model=boats_model.get(boatnumber)
objectif_list=[cap_bonne_esperance_pos,Leuuwin_pos,horn,canary,portugal,goal]
bestroad=[]
if nb_dev>5:
print("maximum 5 deviation")
return None
else:
for goal_number in range(next_point, len(objectif_list)):
if goal_number>3:
spread_size=0.6
else:
spread_size=4
road,reach_goal=from_a_to_b_with_random_variation(start_post,objectif_list[goal_number],cy_Uwind,cy_Vwind,ly_Uwind,ly_Vwind,start_date,boat_speed_model,nb_random_mid_point,spread_size=spread_size, nb_dev=nb_dev, first=True)
if not reach_goal:
road,_=from_a_to_b_with_random_variation(start_post,objectif_list[goal_number],cy_Uwind,cy_Vwind,ly_Uwind,ly_Vwind,start_date,boat_speed_model,nb_random_mid_point,spread_size=spread_size, nb_dev=nb_dev, first=True,use_deg_special_untrapper=True)
bestroad+=road
start_post=bestroad[-1]
#print(start_post)
return bestroad
def from_a_to_b_with_random_variation(a, b, Unc_file_cy, Vnc_file_cy, Unc_file_ly, Vnc_file_ly, start_date,
boat_speed_model, nb_random,spread_size=1, nb_dev=1, first=False, max_len=None,
use_deg_special_untrapper=False):
#function use to search the path between two balise
bestroad, not_trapped=from_a_to_b(a,b,cy_Uwind,cy_Vwind,ly_Uwind,ly_Vwind,start_date,boat_speed_model,spred_want=spread_size, nb_dev=nb_dev,use_deg_special_untrapper=use_deg_special_untrapper)
changed=False
reach_goal=True
if nb_random>0 and len(bestroad)-1>6:
random_place_on_path=int(np.round(len(bestroad)/2,0))
randoms_mid_points=get_ecart_to_line(a,b,nb_dev)
for new_point in randoms_mid_points:
#print(new_point)
bestroad_to_mid, _=from_a_to_b_with_random_variation(a,new_point,cy_Uwind,cy_Vwind,ly_Uwind,ly_Vwind,start_date,
boat_speed_model, nb_random-1,spread_size=spread_size, nb_dev=nb_dev, first=False,
max_len=len(bestroad),use_deg_special_untrapper=use_deg_special_untrapper)
bestroad_from_mid, reach_goal_from=from_a_to_b_with_random_variation(bestroad_to_mid[-1],b,cy_Uwind,cy_Vwind,ly_Uwind,ly_Vwind,
start_date+len(bestroad_to_mid)*6,boat_speed_model, nb_random-1,
spread_size=spread_size, nb_dev=nb_dev,first=False, max_len=len(bestroad),use_deg_special_untrapper=use_deg_special_untrapper)
new_road=bestroad_to_mid+bestroad_from_mid
if (len(new_road)<=len(bestroad) and reach_goal_from) or ((not not_trapped) and reach_goal_from):
bestroad=new_road
changed=True
if (not not_trapped) and (not changed):
br=[]
for item in bestroad:
if item not in br:
br.append(item)
bestroad=copy.deepcopy(br)
reach_goal=False
return bestroad, reach_goal
def transform_a_b_to_same_base(a,b):
#This function is used when going from a positive longitude to a negative one
new_a=copy.deepcopy(a)
new_b=copy.deepcopy(b)
if abs(a[1]- b[1])>180:
if a[1]<0:
new_a=(new_a[0],360+a[1])
else:
new_b=(new_b[0],360+b[1])
return new_a, new_b
def from_a_to_b(a, b, Unc_file_cy, Vnc_file_cy, Unc_file_ly, Vnc_file_ly, start_date,boat_speed_model, spred_want=6, nb_dev=1, max_len=None,use_deg_special_untrapper=False):
#search the best way to go from a to b using only the current wind data to search which option is the best at each date
path=[]
path.append(a)
#get lon, lat, and time in file. Reindex the ones in the last year to make them be used while don't have anymore data
lon_in_file,lat_in_file, time_in_file_cy=get_lon_lat_time_from_nc(cy_Vwind)
time_in_file_cy-=time_in_file_cy[0]
_, _, time_in_file_ly=get_lon_lat_time_from_nc(ly_Uwind)
time_in_file_ly-=time_in_file_ly[0]
for t in range(len(time_in_file_cy)):
time_in_file_ly[t]+=max(time_in_file_cy)
if max_len==None:
max_len=1000
quarter_day=0
not_trapped=True
while get_distance(a,b)>spred_want and quarter_day<max_len and not_trapped:
windspeed, windangle=get_values_wind(Unc_file_ly, Vnc_file_ly, lat_in_file, lon_in_file, time_in_file_ly, a[0], a[1], start_date )
options=get_options(a,b,windangle, windspeed, boat_speed_model, nb_dev=nb_dev,use_deg_special_untrapper=use_deg_special_untrapper)
dist_tmp=[]
new_points=[]
for opt in options:
endpoint=where_i_arrive(a, opt, b=b)
if endpoint[0]>-90 and endpoint[0]<90:
if check_if_line_on_ocean(a, endpoint):
new_points.append(endpoint)
esb, bsb=transform_a_b_to_same_base(endpoint, b)
dist_tmp.append(get_distance(esb, bsb))
if len(new_points)>0:
a=new_points[np.argmin(dist_tmp)]
path.append(a)
start_date+=6
quarter_day+=1
else:
not_trapped=False
path+=[a for iter in range(1000)]
asb,bsb=transform_a_b_to_same_base(a,b)
if get_distance(asb,bsb)>spred_want:
not_trapped=False
return path, not_trapped
def get_options(a,b,windangle, windspeed, boat_speed_model, nb_dev=1,use_deg_special_untrapper=False):
#get the options the boat have at each date
a, b = transform_a_b_to_same_base(a, b)
base_vect=np.array((b[0]-a[0],b[1]-a[1]))/get_distance(a, b)
val=np.arccos(base_vect[0])
dev=np.pi/2/(nb_dev+1)
if use_deg_special_untrapper:
dev=np.pi/(nb_dev+1)
options=[]
for dev_i in np.arange(-nb_dev, nb_dev+1,1):
deg=get_angle_wind_boat(windangle,math.degrees(val+dev*dev_i))
#print("speed: "+str(windspeed)+" and deg: "+str(deg))
boat_speed=boat_speed_model.predict(poly.fit_transform([transform_data_for_model(windspeed, deg)[0]])).reshape(-1,1)
boat_distance_done=boat_speed[0][0]/10
options.append((np.cos(val+dev*dev_i)*boat_distance_done,np.sin(val+dev*dev_i)*boat_distance_done))
return options
def check_if_line_on_ocean(a, b):
is_on_ocean=True
for iter in range(0,7):
c=((a[0]*iter+b[0]*(6-iter))/6,(a[1]*iter+b[1]*(6-iter))/6)
if not globe.is_ocean(c[0], c[1]):
is_on_ocean=False
return is_on_ocean
def get_ecart_to_line(a_r, b_r, nb_ecart):
#create deviation from the most straight road to force the function to test other way of reaching point
if np.sign(b_r[1])!=np.sign(a_r[1]):
if a_r[1]>0:
a,b= a_r, (b_r[0],360+b_r[1])
else:
a,b=(a_r[0],360+a_r[1]), b_r
transformed=True
else:
a,b=a_r, b_r
transformed=False
points=[]
vect_a_b=(b[0]-a[0],b[1]-a[1])
D=(a[0]+vect_a_b[0]/2,a[1]+vect_a_b[1]/2)
AD_dist=get_distance(a,D)
deg=np.sin(np.pi/4)
C=((a[0]+b[0])/2,(a[1]+b[1])/2)
if transformed and C[1]>180:
C=(C[0],C[1]-360)
points.append(C)
for deg_num in range(nb_ecart):
deg=np.pi/3/7*(deg_num+1)
DC_dist=AD_dist*np.sin(deg)/np.sqrt(1-np.sin(deg)**2)
DC_vect=(vect_a_b[1]*DC_dist/AD_dist,vect_a_b[0]/2*DC_dist/AD_dist)
C=(D[0]+DC_vect[0], D[1]+DC_vect[1])
if transformed and C[1]>180:
C=(C[0],C[1]-360)
if C[0]>-90 and C[0]<90 and C[1]>-180 and C[1]<-180:
if globe.is_ocean(C[0],C[1]):
points.append(C)
else:
deg=np.pi/3/7/14*(deg_num+1)
DC_dist=AD_dist*np.sin(deg)/np.sqrt(1-np.sin(deg)**2)
DC_vect=(vect_a_b[1]*DC_dist/AD_dist,vect_a_b[0]/2*DC_dist/AD_dist)
C=(D[0]+DC_vect[0], D[1]+DC_vect[1])
if transformed and C[1]>180:
C=(C[0],C[1]-360)
if C[0]>-90 and C[0]<90 and C[1]>-180 and C[1]<-180:
if globe.is_ocean(C[0],C[1]):
points.append(C)
deg=-np.pi/3/7*(deg_num+1)
DC_dist=AD_dist*np.sin(deg)/np.sqrt(1-np.sin(deg)**2)
DC_vect=(vect_a_b[1]*DC_dist/AD_dist,vect_a_b[0]/2*DC_dist/AD_dist)
C=(D[0]+DC_vect[0], D[1]+DC_vect[1])
if transformed and C[1]>180:
C=(C[0],C[1]-360)
if C[0]>-90 and C[0]<90 and C[1]>-180 and C[1]<-180:
if globe.is_ocean(C[0],C[1]):
points.append(C)
else:
deg=-np.pi/3/14*(deg_num+1)
DC_dist=AD_dist*np.sin(deg)/np.sqrt(1-np.sin(deg)**2)
DC_vect=(vect_a_b[1]*DC_dist/AD_dist,vect_a_b[0]/2*DC_dist/AD_dist)
C=(D[0]+DC_vect[0], D[1]+DC_vect[1])
if transformed and C[1]>180:
C=(C[0],C[1]-360)
if C[0]>-90 and C[0]<90 and C[1]>-180 and C[1]<-180:
if globe.is_ocean(C[0],C[1]):
points.append(C)
return points
def where_i_arriveold(start, move, b=None):
#function to calculate the reach point from a start point and a move calculated before with wind and boat datas
if b==None:
new_lon=start[1]+move[1]
if new_lon>=180:
new_lon-=360
elif new_lon<=-180:
new_lon+=360
reach_point=tuple([start[0]+move[0], new_lon])
else:
if b[1]-start[1]<move[1] and b[0]-start[0]<move[0]:
reach_point=b
else:
new_lon=start[1]+move[1]
if new_lon>=180:
new_lon-=360
elif new_lon<=-180:
new_lon+=360
reach_point=tuple([start[0]+move[0], new_lon])
return reach_point
def where_i_arrive(start, move, b=None):
#function to calculate the reach point from a start point and a move calculated before with wind and boat datas
new_lon=start[1]+move[1]
if new_lon>=180:
new_lon-=360
elif new_lon<=-180:
new_lon+=360
reach_point=tuple([start[0]+move[0], new_lon])
if b!=None:
esb, bsb=transform_a_b_to_same_base(start, b)
esn, bsn=transform_a_b_to_same_base(start, reach_point)
if get_distance(esb, bsb)<get_distance(esn, bsn):
reach_point=b
return reach_point
def get_distance(a, b):
#return the euclidian distance between two points
return np.sqrt((b[0]-a[0])**2+(b[1]-a[1])**2)
def get_angle_wind_boat(windangle, degree):
angle=abs(windangle-degree)
if angle>180:
angle-=90
if angle>180:
angle=angle-180
return(abs(angle-180))
Désormais nous pouvons pour chaque bateau calculé un chemin théorique pour atteindre l'arrivé en utilisant ces fonctions
####Choose parameters here ### Caution, high value mean very long calculation time ###
nb_dev=5
nb_ecart=5
#############################################################################
def modify_values_chemin_for_graph(chemin):
pass_max_lon=False
for pos in range(1,len(chemin)):
if not pass_max_lon and chemin[pos-1][1]>0 and chemin[pos][1]<0:
pass_max_lon=True
if pass_max_lon:
chemin[pos]=(chemin[pos][0],chemin[pos][1]+360)
return chemin[1:]
if launch_long_calculations:
df_with_new_values=copy.deepcopy(df_complete)
df_with_new_values=df_with_new_values[['Voile', 'Skipper', 'position','date']]
df_with_new_values['date']=df_with_new_values['date'].map(transform_date)
start_date=np.sort(df_with_new_values['date'].unique())[-1]
for boat in boats_numbers:
tmp=df_with_new_values[df_with_new_values['date']==start_date]
start_pos=tmp[tmp['Voile']==boat]['position'].values[0]
skipper=tmp[tmp['Voile']==boat]['Skipper'].values[0]
chemin=find_best_way_to_goal(cy_Uwind,cy_Vwind,ly_Uwind,ly_Vwind,start_date,boats_model_lasso, boat,
nb_ecart, start_pos, next_point=0,nb_dev=3)
chemin=modify_values_chemin_for_graph(chemin)
chemin_date=[start_date+6*iter for iter in range(1,len(chemin)+1)]
Voile_current_chemin=[boat for iter in range(len(chemin))]
Skipper=[skipper for iter in range(len(chemin))]
df_tmp=pd.DataFrame({'Voile': Voile_current_chemin, 'Skipper': skipper, 'position': chemin, 'date': chemin_date})
df_with_new_values=df_with_new_values.append(df_tmp)
delete_if_exist_and_recreate_directory(os.getcwd()+"\\temporary_save_file\\")
df_with_new_values.to_csv(os.getcwd()+"\\temporary_save_file\\result_to_end_run.csv")
else:
df_with_new_values=pd.read_csv(os.getcwd()+"\\temporary_save_file\\result_to_end_run.csv")
df_with_new_values['position']=df_with_new_values['position'].map(lambda x: eval(x))
Et sortir un classement final
tmp=df_with_new_values.groupby(['Skipper']).max().sort_values(by='date')
tmp=tmp.join(df_complete.groupby(['Skipper'])['Foil'].unique().map(lambda x: x[0]))
tmp=tmp.join(df_complete.groupby(['Skipper'])['Année lancement'].unique().map(lambda x: x[0]))
tmp=tmp.join(df_complete[df_complete['date']==df_complete.groupby(['Skipper']).max()['date'].values[0]].groupby(['Skipper']).max()['Rang'])
tmp=tmp.rename(columns={'Rang': 'Dernier Classement'})
tmp['durée parcours (jours)']=tmp['date'].map(lambda x: (x-7502)/24)
tmp=tmp.reset_index().reset_index().rename(columns={"index": "Rang simulé"})[['Skipper','Voile','Rang simulé','Dernier Classement','Foil','Année lancement','durée parcours (jours)']]
tmp['Rang simulé']=tmp['Rang simulé'].map(lambda x: x+1)
tmp
On voit que le modèle est quelque peu optimiste avec près de 20 concurent dépassant le record du Vendée Globe ! Cela est explicable par le fait qu'il ne prend pas en compte les autres conditions météorologique, et que la vitesse des bateaux est croissante de la vitesse des vents, il cherchera donc a trouver les zones de vents les plus rapide. Dans la réalité les skippers feront parfois le choix de ne pas prendre les zones les plus venteuse du fait des dangers que cela peut représenter pour le navire et eux même. Enfin il ne prend pas en compte les avaris possible, mais surtout utilise des données de vent fixé à l'avance. Dans la réalité les prévisions s'étalent au maximum sur les 16 prochains jours et il est donc difficiles de comparer une route tracer dans ces conditions à celle de notre algorithme. Il sera toutefois interressant de regarder les trajets emprunté dans la partie suivante.
Partie Visualisation:¶
Afin de pouvoir visualiser au mieux la course j'ai choisi d'utiliser un graphe animé à l'aide de plotly express. Ce package permet a la fois d'inserer simplement les données sur la carte, mais permet aussi d'obtenir un graphe interactif contenant les informations des bateaux, et d'animer celui-ci.
#Prepare a dataframe for the graph
df_for_graph_tmp=copy.deepcopy(df_complete)
df_for_graph_tmp['lon']=df_for_graph_tmp['position'].map(lambda x: x[1])
df_for_graph_tmp['lat']=df_for_graph_tmp['position'].map(lambda x: x[0])
dates=np.sort(df_for_graph_tmp['date'].unique())
df_for_graph=df_for_graph_tmp[df_for_graph_tmp['date']==dates[0]]
df_for_graph['time_index']=0
for iter in range(1,len(dates)):
dt=dates[:iter]
df_tmp=df_for_graph_tmp[df_for_graph_tmp['date'].map(lambda x: True if x in dt else False)]
df_tmp['time_index']=iter/4
df_for_graph=df_for_graph.append(df_tmp)
Creation du graphe sur les données réelles¶
#Create the graph, get patient can take up to one minute because of the animation
fig = px.line_mapbox(df_for_graph,
lat='lat',
lon='lon',
color="Skipper",
animation_frame = 'time_index',
hover_data=['time_index','Voile'],
zoom=1,
height=800,
title = 'Visualization of the Vendee Globe par Skypper')
fig.update_layout(
mapbox_style="white-bg",
mapbox_layers=[
{
"below": 'traces',
"sourcetype": "raster",
"sourceattribution": "United States Geological Survey",
"source": [
"https://basemap.nationalmap.gov/arcgis/rest/services/USGSImageryOnly/MapServer/tile/{z}/{y}/{x}"
]
}
],
mapbox_center_lat=0,
margin={"r":0,"t":0,"l":0,"b":0},
)
fig.update_layout()
fig.show()